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Abstract. We present a new method for determining Weakly Interacting Massive Particle 
(WIMP) properties in future tonne scale direct detection experiments which accounts for 
uncertainties in the Milky Way (MW) smooth dark matter distribution. Using synthetic data on 
the kinematics of MW halo stars matching present samples from the Sloan Digital Sky Survey, 
complemented by local escape velocity constraints, we demonstrate that the local dark matter 
density can be constrained to ~ 20% accuracy. For low mass WIMPs, we find that a factor 
of two error in the assumed local dark matter density leads to a severely biased reconstruction 
of the WIMP spin-independent cross section that is incorrect at the 15cr level. We show that 
this bias may be overcome by marginalizing over parameters that describe the MW potential, 
and use this formalism to project the accuracy attainable on WIMP properties in future 1 ton 
Xenon detectors. Our method can be readily applied to different detector technologies and 
extended to more detailed MW halo models. 
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1. Introduction 

The detection of dark matter via elastic scatterings in underground detectors is a foremost goal 
of experimental physics in the coming years. Upper limits on the spin-independent WIMP- 
nucleon cross section E l2l El are now beginning to carve into the parameter space of the 
constrained Minimal Supersymmetric Standard Model (cMSSM) dH [Sj, which provides a 
well-motivated theoretical framework for the dark matter. Favored MSSSM parameter space 
is expected to be fully probed as future detectors reach the tonne scale. 

A robust interpretation of the experimental limits (or signals) from detectors requires an 
understanding of uncertainties associated with the distribution of dark matter in the Milky 
Way (MW), specifically in the vicinity of the solar neighborhood. Though long regarded a 
nuisance in the prediction of scattering event rates [6], the potential for the uncertainty in the 
local dark matter distribution to bias the reconstruction of the WIMP mass and cross section 
has yet to be fully quantified. 

Here we study how well, in a realistic scenario, forthcoming direct detection experiments 
can do in reconstructing the WIMP mass and cross section, accounting for an uncertain 
smooth component of the dark matter distribution, which is independently constrained using 
spectroscopic observations of distant halo stars [j7ll8| and measurements of the escape speed 
near the solar circle BU. We quantify the bias introduced in the reconstruction of WIMP 
properties when the uncertainty in the MW halo model is neglected. The results we present are 
specifically focused on the spin-independent WIMP-nucleon cross section, but are otherwise 
independent of the specific particle physics framework. 

2. Milky Way Halo 

Rates in direct detection experiments are proportional to the dark matter density at the solar 
radius. However the potential at the solar radius is likely dominated by baryons, so a self- 
consistent approach must model the combined potentials of baryons and dark matter. Disk 
stars on circular orbits trace the local potential fTOl, and the dark matter contribution may be 
determined via continuation of constraints at large radii and under the assumption of a smooth 
dark matter density profile. This approach is of course subject to systematic uncertainties in 
the parameterization of the halo potential, though it provides a starting point for understanding 
the relative contributions from each of the key MW mass components. As we argue here, even 
a simplified MW model already represents a very considerable improvement in bringing under 
control systematic errors in reconstructed WIMP properties. 

2.1. Model for Mass Distribution 

We consider a MW halo model that includes a central bulge, disk, and dark halo ©[mi- The 
bulge is modeled as a spherically- symmetric Hernquist potential, 0buige('') = -GMbuige/('' + 
Co), with Co ~ 0.6 kpc and a total mass Mbuige = l.SxlO^^M©. The disk surface density is taken 
of the form S(i?) oc e ^^/^**^ with R the cylindrical coordinate for the axially- symmetric disk 
and its scalelength is Z^disk - 4 kpc, which is included as a free parameter in our conservative 
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model (see below). The potential from this component is clearly non- spherical, though it 
may be fairly accurately modeled by a spherical distribution that has the same mass interior 
to a Galactocentric radius r = R. The spherically- symmetric disk potential is taken as 
0disk = -GMdisk(l - e"''^'''"*)/r, where the total disk mass is Mdisk = 5 x 10"' Mq. The peak 
circular velocity for the spherical fit to the disk potential at ~ 2Z?disk is ~ 15% less than the 
circular velocity of the exponential disk and asymptotes to < 5% of the true mass distribution 
for large radii [fT2ll . For the dark matter halo we take a five parameter model, 

^^'^ = (r/ro)''[l+^r/ro)^](-«)/^- 
The escape velocity, Vesc, and the circular velocity at the solar radius, vq, are determined from 
the combined potential of the three components. 

As is the case with any parametric model, the results we present will likely vary if the 
model is not an accurate description of the true Galaxy. As such, we view our analysis as 
means of estimating uncertainties and bias on key dark matter parameters for a well-defined, 
though perhaps simplified, MW model, and acknowledge that the above model provides a first 
step towards consideration of a wider range of Galactic models in this context. For example, 
to better model observations lfT3]| triaxial models for the MW halo may be considered which 
increase the dark matter density in the disk. However, on the theoretical side there are 
uncertainties to the predicted MW halo shape, e.g. simulations suggest that adding gas cooling 
tends to make halos more spherical [[T4l|. 



2.2. Stellar Kinematics 

From the above parametric model for the Galaxy, we simulate stellar kinematics data which 
we then employ to determine the accuracy with which the particle properties of the dark matter 
and its astrophysical distribution can be reconstructed. Our synthetic data sample consists of 
2000 stars distributed uniformly with Galactocentric radii in the range 5-40 kpc, similar to 
the distribution in recent Sloan Digital Sky Survey measurements IfTTll . The distance, d, from 
a halo star to the Sun is given by = -I- r| - 2 r cos 6, so that the line-of-sight velocity 
to the star is ^i' = r(r/d - cos 9/d) + rrQO sin 6/d. Here 6 is the spherical polar angle and 
= 8.5 kpc is the distance from the Sun to the Galactic center, which we keep fixed. The 
line-of-sight velocity dispersion is defined as the average of d^ over the solid angle in Galactic 
coordinates, cr^^^ = (d^). Performing this average gives [il5il . 



(Tios = ^Jl-A{r)/3{r)cr,, (2) 

where 

r + ro 



2,2 / 2 _ 2^2 



r-rQ 

.2/^2 



(3) 



The velocity anisotropy parameter is J3{r) = 1 - o-g/crf, which we assume to be constant as a 
function of radius. 

Assuming that the dark matter halo and the tracer population of halo stars are in 
equilibrium, the radial jeans equation for the intrinsic, or statistical, dispersion cr,. is 
d(p^o-j) Ip^ajfiir) d(p 

— + = -Gp*—, (4) 

dr r dr 
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where p*(r) is the spatial density profile for the tracer population. The potential on the 
right-hand side is the sum of the potential from the disk, bulge, and the dark matter, 
= 0disk + 0biuge + 0dm- The density of halo stars is modeled to fall off" according to p^, oc r~^-^ 
beyond 10 kpc and remain constant for r < 10 kpc [il6il . 

3. Direct Detection 

To calculate WIMP-nucleon recoil event rates, we must estimate both the local WIMP density 
and the velocity distribution of the dark matter. For the latter, we appeal to recent results from 
high resolution numerical simulations, which find only mild deviations from an isotropic and 
Maxwellian velocity distribution in the solar neighborhood IfTTl [T8l [T9l . This implies 10% 
deviations in recoil event rates relative to the standard isotropic Maxwellian assumption [|20l . 
Motivated by these results, we assume an isotropic Maxwellian form for the dark matter 
velocity distribution, characterized by Vesc, vo, and local density p(8.5kpc) f2T\, which are all 
computed from our MW model described above and constrained as explained below. 

We consider the case of a Xe target, and assume everywhere a 1 tonne-yr exposure. In 
order to model the energy resolution of the detector, the recoil energy spectrum is smoothed 
with an energy-dependent Gaussian window function of standard deviation asiE) = 1.5 Ve. 
For definiteness, our results below are obtained by considering 5 equally- spaced energy bins 
within the recoil range 5-30 keV, and adopting an independent Poisson likelihood in each 
bin, with negligible background counts. However, our results are insensitive to the details 
of the binning scheme and to the precise choice of the upper cutoff" in the recoil spectrum. 
The nuclear form factor for spin-independent scatterings is taken to be the standard Helm 
form [|2T|. For the sake of simplicity, we do not include the yearly variation in the event rate 
resulting from the projection of the Earth's orbital velocity on the Galactic plane; we find that 
including this eff'ect has negligible eff'ect on the results we present. 

4. Statistical methodology 

We consider the following fiducial parameters for the Milky Way: = [po, ro,/3, a, b, c, ^disk) 
= {7.4 X 10^ Mq kpc"^, 20 kpc, 0, 1, 1, 3, 4 kpc}. This set of dark halo parameters is broadly 
motivated by numerical simulations of Milky Way mass halos [|22ll . and as discussed above 
the mass and scalelength of the disk, Z^disk, are of the order observed for the Milky Way [|23l . 
For this set of parameters the local dark matter density is 0.32 GeV cm"^, the escape velocity 
Vesc = 547 km s"^ and the circular velocity vq = 209 km s"^ . We find that the results presented 
below are insensitive to similar sets of fiducial model parameters, provided that they roughly 
match both the local dark matter density and the mass and circular velocity of the Milky Way. 
The likelihood is centered around the fiducial values, and is the product of three components: 
the projected direct detection constraints, the projected stellar kinematics constraints, and 
present-day escape velocity measurements. 

We define the likelihood for the projected velocity component along the line-of-sight 
to a given star to be Gaussian with zero mean and a variance given by the ctios. The total 
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likelihood is then obtained by multiplication of the Gaussian likelihood for each star. The 
dispersion in the likelihood is dominated by the intrinsic dispersion at a given radius, i.e. we 
ignore contributions to the likelihood that involve uncertainty on the measured velocity of 
each star. Assuming that the model parameters only enter through the velocity dispersion, it 
is straightforward to derive the Fisher information matrix, which is defined as the second 
derivative of the log of the likelihood function, with respect to the model parameters, 
(d^JL/d&ad&b)- For our likelihood the Fisher matrix is simply given by 

1=1 loS,! " " 

Written in this form we assume that the mean velocity does not depend on our model 
parameters; this is likely to be a good assumption if we consider that the intrinsic dispersion 
of the halo stars dominates the contribution due to ordered rotational motion. 

We take the constraint on the escape velocity to be Gaussian with mean Vgsc = 544 km 
s'^ and standard deviation cr = 33 km s~' |l9l. For simplicity we do not consider additional 
constraints on the MW mass distribution that come from the inner and outer rotation curves or 
from the motions of stars in the solar neighborhood, as considered in Refs. [f24ll23ll25l . Our 
motivation for including only the halo stars and the escape velocity constraints in our data set 
is that these quantities are the most direct tracers of the MW dark matter halo. Given these 
inputs, the projected constraints on both the local dark matter density and the reconstructed of 
WIMP properties we present below should be viewed as a conservative estimate of the errors 
attainable on these parameters. 

We adopt a Metropolis-Hastings algorithm to perform a Markov chain Monte Carlo scan 
over the parameters defining the MW halo model described above, augmented by the WIMP 
mass and its spin-independent scattering cross section cr^^. We assume fiat priors on all 
MW halo parameters over the ranges {log[po/(Mokpc~^)], ro/kpc,jS, a, b, c, ^disk/kpc} = {[5 : 
8], [10 : 60], [-1 : 1], [0 : 2], [0 : 2], [2 : 4], [2 : 6]}, as well as flat priors over sufliciently 
wide ranges of log m^^ and log cr^p'. We find that our results are insensitive to the details of the 
prior choice. 

5. Results 

We study three representative cases, which differ in the number of halo parameters that are 
marginalized over. Om fixed model takes to be fixed to an assumed value, which might 
not match the correct one. This corresponds to neglecting all astrophysical uncertainties 
while not including any of the astrophysical constraints. Our baseline model varies the 
following three parameters: {logpo, rQ,/3}, which are constrained using the astrophysical data 
and marginalized over when considering direct detection constraints, while the remaining 
four parameters are fixed to their fiducial value. In addition to the three parameters above, our 
conservative model also marginalizes over the remaining four parameters of the halo model, 
{a, b, c, bdisk). 

Figure \T\ shows the ID posterior marginal probabilities for several parameters for the 
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Figure 1. Reconstruction of 3 halo model parameters (logpo, '"o,y6) and constraints on WIMP 
scattering rates-related quantities at the solar circle (p(8.5 kpc), Vesc, t^'o) using 2000 halo stars 
and a local determination of the escape velocity. Green (tighter) and blue (wider) curves show 
the results for the 3 parameter (baseline) and 7 parameter (conservative) halo model. The red 
diamond is the true value. 



baseline and the conservative models, after the halo stars and escape velocity constraints 
have been applied. The distributions are correctly peaked around the true value and in 
particular the density of dark matter at 8.5 kpc can be constrained to 22% (37%) at Icr in the 
baseline (conservative) model. We find that the local escape velocity constraint is important 
in tightening the bounds on the WIMP local density. 

Figure |2] depicts the reconstruction of and cr^' for a 1 tonne-yr Xe exposure for 
the baseline (filled/green) and conservative (filled/blue) models, assuming a fiducial value 
CTp' = 10"^ pb. We show the case of a low mass (m^ = 50 GeV, left panel) and a high 
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Figure 2. Reconstruction of m^, cr^,' under various assumptions: dark matter halo parameters 
fixed to assumed values (solid, black), marginalizing over baseline halo model (filled/green), 
marginalizing over conservative halo model (filled/blue). In all cases inner and outer contours 
represent 68% and 95% c.l. limits. The red diamond gives the true value. The left panel is for 
a 50 GeV WIMP mass, the right panel assumes a 500 GeV WIMP (In the right panel, we only 
show 68% c.l. for the case of fixed galactic parameters for clarity). The lower (upper) solid 
black contours illustrate the bias in the reconstruction assuming incorrect values for the local 
dark matter density a factor of 2 above (below) the true value. 



mass WIMP (m^ = 500 GeV, right panel), to highlight the difference between WIMP mass 
scales that are well-constrained and poorly constrained ll26l l27l l28l . At = 50 GeV, the 
degeneracy between m^^ and cr^^ is broken by recoil energy spectral information; we find that 
a minimum of 3 evenly spaced bins suffices to break the degeneracy. The fiducial model with 
= 50 GeV (m^ = 500 GeV) produces 383 (79) events within our energy window for a 
tonne-yr exposure for Xe. 

Figure |2] also shows the bias in the reconstruction resulting from assuming the incorrect 
local dark matter density. For a 50 GeV WIMP, fixing the local dark matter density a factor 
of two above or below the true value biases the reconstruction of the WIMP cross section 
by ~ 15cr. The effect is less severe for a 500 GeV WIMP, of order ~ 5cr, though this case 
is intrinsically much less constrained. This systematic bias would not be detectable by the 
usual goodness-of-fit test, for in all three cases depicted in Figure |2] the reduced chi-square 
of the best-fit point is statistically indistinguishable. Figure |2] clearly shows the necessity of 
accounting for uncertain astrophysics in the determination of /n^, cr^', even in a simplified 
MW halo model. For the well constrained 50 GeV WIMP case, marginalizing over the halo 
model parameters correctly recovers the true point with no bias, while the error in the mass 
increases from 9% to 14% (24%) for the baseline (conservative) model. The impact on 
the accuracy on log cr^^ is more severe, increasing the relative error by a factor of 5 under 
the baseline model and by a factor of 20 for the conservative case. However, this reduced 
statistical accuracy is compensated by a greatly increased robustness to systematic bias. 



Reconstructing WIMP Properties in Direct Detection Experiments 



8 



6. Conclusion 

We have presented an analysis of the constraints attainable on WIMP parameters using future 
direct detection experiments. Our key results are twofold: first we have quantified how 
well WIMP properties can be constrained for an optimistic particle mass in tonne scale 
detectors, and second we have demonstrated that including an explicit model of the MW halo 
can overcome a potentially severe (~ 15cr) bias in the reconstruction of WIMP parameters 
when the wrong fiducial values for the halo are assumed. In future analyses our MW halo 
model can be refined to include triaxial shapes for the dark and luminous components. The 
contribution from dark matter substructure If29l [TSlI or a rotating component BOi may be 
considered. Further one may account for the non-Maxwellian velocity distribution OTll . and 
a multi-component spectral fit to the WIMP and astrophysical background spectra may be 
incorporated [l32l[331 . An analysis along these lines will be crucial to interpret the limits and 
measurements from forthcoming direct detection experiments. 
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